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ABSTRACT 


A least-squares  power-spectrum  analysis  of  122  years  of  Zurich  dally 
sunspot  numbers  yields  a statistically  significant  peak  at  12.07151.002 
days  period.  This  feature  of  the  sunspot  spectrum  may  be  associated 
with  the  peak  at  12.22  days  (sideral)  which  Dicke  (1976)  found  in  his 
oblateness  data,  and  may  be  attributable  to  the  sun's  core  if  it  rotates 
at  either  12.0715  days  or  24.1430  days  period  (synodic). 


I.  Introduction 

The  differential  rotation  of  the  solar  surface  is  well  established, 
but  different  measurement  techniques  give  different  rates  (Gilman  1974, 
Howard  1976)  and  different  patterns,  some  "tracers"  showing  rigid  rota- 
tion or  very  little  differential  rotation  (Uilcox  et  al . 1970,  Timothy 
et  al . 1975,  Adams  1976).  There  is  also  some  evidence  that  the  surface 
rotation  measured  by  the  same  technique  is  not  constant  in  time  (Howard 
1976).  There  is  at  the  present  time  no  single  widely  accepted  theoreti- 
cal explanation  of  solar  differential  rotation. 

If  the  rotation  of  the  sun's  surface  is  poorly  understood,  it  is  not 
surprising  that  we  know  even  less  about  its  rotation  below  the  visible 
layers.  Various  hypothesis  about  the  rate  and  character  of  the  internal 
rotation  of  the  sun  (see  for  example  Oicke  1964,  Schatten  1977)  have 
been  proposed,  and  sophisticated  observational  techniques  and  data  ana- 
lysis may  soon  yield  information  about  the  rotation  of  the  solar  inte- 
rior (Deubner  et  al . 1978).  It  seems  likely  that  the  differential  rota- 
tion is  confined  to  the  convective  region  of  the  sun,  but  the  extent  of 
this  region  is  unknown.  It  seems  not  unreasonable  to  expect  rigid  rota- 
tion of  the  sun's  radiative  core  and  perhaps  also  some  inner  portion  of 
its  convection  zone. 

Dicke  has  reported  evidence  for  some  rotation  more  rapid  than  the 
usual  surface  rotation  rates  (Dicke  1974,  1976).  Dicke's  most  recent 
view  (1976)  seems  to  be  that  there  exists  a photospheric  perturbation 
rotating  with  a 12.22  day  (sidereal)  period.  As  an  independent  check  of 
the  reality  of  the  rotation  reported  by  Dicke  (1974,  1976),  we  have  gen- 
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» crated  least-squares  spectra  of  44520  daily  sunspot  numbers  and  find  a 

peak  at  12.0715±0.002  days,  near  Dicke's  suggested  rotation  period  of 
12.64  days  (synodic).  The  statistical  significance  of  this  peek  and  its 
relation  to  Dicke's  results  are  discussed  in  Section  2.  The  results  of 
the  analysis  of  Section  2 are  summarized  and  some  possible  interpreta- 
tions are  briefly  discussed  in  Section  3. 
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II.  Spectrum  Analysis 


He  have  analyzed  the  run  of  Zurich  daily  sunspot  numbers  (see  Mald- 

meier  1961)  extending  from  January  7,  1649  to  November  26,  1970  by  a 

least-squares  procedure.  Me  minimize 

^ 2 

V3  = - a - b cos - n)  - c sin(-^  n)  , (2.1) 

n=l  L J 


where  Xn  is  the  sunspot  number  and  n counts  days  from  1 to  44520,  vary- 
ing a,  b and  c for  a range  of  values  of  the  period  P.  Then  S = b2+c2 
provides  an  estimate  of  the  spectral  power  at  period  P.  The  sums  re- 
quired to  calculate  the  least-square  fits  were  generated  using  a 131072 
element  FFT . Figure  1 displays  the  estimated  spectral  power  as  a func- 
tion of  period  for  periods  greater  than  4 days.  The  large  number  of 
spectral  estimates  (32768)  precludes  plotting  all  the  individual  points 
so  the  spectrum  has  been  grouped  into  64  frequency  intervals  each  con- 
taining 512  spectral  estimates.  The  median,  70th,  90th  and  99th  percen- 
tiles for  each  of  the  64  bins  are  displayed  in  Figure  1.  The  sun  symbol 
in  Figure  1 shows  the  peak  at  12.0715  days.  There  is  obviously  a great 
deal  of  power  in  the  very  low  frequency  portion  (associated  with  the 
eleven  year  solar  cycle)  and  in  the  range  of  the  familiar  surface  rota- 
tion periods. 

Me  wish  to  assess  the  significance  of  the  peak  at  12.0715  days.  One 
common  method  used  to  assess  the  significance  of  a least-squares  fit  is 
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the  F test  (see  for  example  Rao  1973).  Me  form  the  statistic 


V -V 

N-3  V1  3 


2,N-3 


(2.2) 


where  Vs  is  defined  by  (2.1)  and  Vi  is  the  sum  of  the  squared  deviations 
from  the  mean.  Me  expect  this  statistic  to  be  distributed  approximately 
as  Fz, 11517  if  the  following  assumptions  are  reasonably  well  satisfied: 
the  residuals,  used  to  form  the  sums  Vi  and  Vs,  (a)  are  distributed  nor- 
mally with  mean  zero  and  the  same  variance,  and  (b)  are  uncorrelated  day 
to  day. 

However,  one  can  see  that  the  sunspot  numbers  are  highly  correlated 
from  day  to  day  by  simple  inspection  of  the  data  (see  Maldmeier  1961). 
Furthermore,  the  daily  sunspot  numbers  are  very  "skewed"  compared  to  a 
sample  selected  from  a normal  population  with  the  same  mean.  It  is  true 
that  the  F test  is  "robust"  (i.e.  it  is  fairly  insensitive  to  the  as- 
sumptions of  normality  being  exactly  satisfied)  but  the  sunspot  number 
data  are  far  from  normally  distributed  and  we  therefore  cannot  expect 
the  F statistic  constructed  using  the  least-square  fits  to  the  raw  sun- 
spot data  to  be  distributed  as  F2, 11517. 

Me  can  see  that  the  F statistic  constructed  from  the  least-square 
spectral  estimates  is  not  distributed  as  F2, 11517  from  the  F statistics 
for  the  high-frequency  (short-period)  portion  of  the  least-square  spect- 
rum. For  the  24538  highest-f requency  least-square  estimates  (periods 
from  4 to  about  16  days),  none  exceeds  the  1%  significance  level.  Me 
have  in  fact  oversampled  in  frequency  by  a factor  of  about  3 to  assure 
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that  no  strong  features  are  neglected,  so  that  there  are  only  about  8000 
independent  spectral  estimates  in  this  range  of  periods.  The  probabil- 
ity that  none  of  8000  independent  F statistics  would  exceed  the  IX  level 
is  (.99),#0°  « 10'3*-*.  Clearly,  the  usual  F test  does  not  provide  an 
adequate  estimate  of  the  significance  of  a particular  spectral  feature. 

Whether  or  not  the  original  data  are  normally  distributed,  the 
least-square  estimates  of  b and  c of  equation  (2.1)  are  effectively  sums 
of  a large  number  of  products  of  individual  data  and  sines  or  cosines 
for  which  the  expectation  values  are  zero,  so  that  they  may  be  expected 
to  be  distributed  normally  with  mean  zero  if  there  are  no  periodic  "sig- 
nals" in  the  data  at  the  frequency  in  question.  Me  therefore  expect  the 
spectral  estimates  to  be  distributed  as  c2x2  with  two  degrees  of  freedom 
(o2x|  ),  where  o2  is  the  variance  (assumed  equal)  of  b and  c.  The  ex- 
pected value  of  the  c2x2  statistic  is  2o2.  There  is  no  simple  way  to 
estimate  what  value  to  expect  for  o;  indeed,  inspection  of  Figure  1 in- 
dicates it  depends  strongly  on  frequency. 

Me  are  primarily  interested  in  periods  near  Dickers  suggested  (sy- 
nodic) rotation  period  of  12.64  days.  Me  expect  that  the  large  broad 
peak  near  27  days  is  real  and  due  to  the  combination  of  a non-uniform 
distribution  of  sunspots  on  the  surface  of  the  sun  and  the  differential 
rotation  of  the  surface.  Me  therefore  will  restrict  our  attention  to 
the  24538  highest  frequencies  (periods  from  4 to  about  16  days).  The 
obvious  trend  in  this  portion  of  the  spectrum  was  removed  by  fitting  a 
quadratic  in  log(f)  to  the  log  of  the  estimated  spectral  power  averaged 


1 

over  10  bins  approximately  equally  spaced  in  log(f).  Having  placed  the 


high  frequency  spectral  estimates  on  approximately  the  same  footing,  ue 


are  in  a position  to  estimate 
the  probability  density. 


the  variance  needed  to  specify 


The  cumulative  density  function  for  a statistic  (in  this  case  the  ad- 
justed spectral  estimates,  S)  distributed  as  ozx£  is 


C = P(S’<S) 


S 

/ 


-S' /2a 

e dS ' = 


-S/2  a‘ 


j . (2.3) 


He  can  estimate  a2  by  comparing  the  theoretical  cumulative  distribution 
with  an  empirical  cumulative  distribution  constructed  from  the  adjusted 
spectral  estimates.  The  spectral  estimates  are  sorted  and  each  of  the 
estimates  is  assigned  a value  of  C according  to 


C = (i-.5)/24538 


(2.4) 


where  i=1  for  the  smallest  adjusted  spectral  estimate  and  i=24538  for 
the  largest. 


Ue  have  used  three  different  methods  to  estimate  o2  from  the  sorted 
adjusted  spectral  estimates.  First  ue  fit  the  empirical  cumulative  dis- 
tribution by  varying  l/2o2  to  minimize 


24538 


S t s* 


+ ln(l-C . 


(2.5) 


and  obtain  e2  = .5025.  Since  we  expect  -1  S;/ln(1-C|)  to  be  approxi- 
mately equal  to  o2,  we  can  average  this  quantity  and  the  reciprocal  of 


it  to  obtain  the  second  and  third  estimates,  which  are  o2  = .4981  and 
o2  = .4996  respectively.  Combining  these  estimates  for  e2  with  the  ad- 
justed power  at  12.0715  days  (S*13.46),  we  estimate  (a  posteriori)  the 
probability  that  this  large  a peak  would  occur  by  chance  to  be 
1.53  x 10'*,  1.35  x 10*‘  or  1.41  x 10'*  respectively. 

As  we  have  indicated,  we  expect  C - (1-e*^/2(2  ).  A plot  of  1 n ( 1—0 
versus  the  sorted  adjusted  spectral  estimates  therefore  should  be  nearly 
a straight  line.  In  Figure  2 the  theoretical  relation  between  (1-0  and 
S for  o2  = .5  is  the  broken  line  and  the  actual  adjusted  spectral  esti- 
mates are  indicated  by  the  solid  curve.  The  points  corresponding  to  the 
five  largest  adjusted  spectral  estimates  are  indicated  by  the  numbers  1 
through  5.  As  expected,  the  empirical  and  theoretical  curves  are  quite 
close  except  for  the  two  largest  adjusted  spectral  estimates  which  are 
associated  with  the  peak  at  12.0715  days.  The  deviation  of  the  empiri- 
cal curve  from  the  theoretical  line  for  (1-C)  between  -10'2  and  -10'*  is 
primarily  due  to  an  excess  of  large  adjusted  spectral  estimates  associ- 
ated with  harmonics  of  the  broad  peak  near  27  days  period. 

Me  conclude  from  the  above  analysis  that  the  adjusted  spectral  esti- 
mates may  reasonably  be  assumed  to  be  distributed  as  c2*2  with  o * .5 
and  that  the  (a  posteriori ) probability  that  the  peak  at  12.0715  is  due 
to  chance  is  -1.4  x 10'*.  However,  we  must  take  account  of  the  fact 
that  the  period  12.0715  was  not  chosen  a priori . but  inferred  from  the 
data. 


Oicke  (1976)  interpreted  the  results  of  his  analysis  as  evidence  for 
some  solar  rotation  with  n period  of  12.641.12  days  (synodic).  Dickers 
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quoted  error  estimates  for  the  period  are  calculated  from  a maximum 


r F 
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likelihood  treatment  of  the  residuals  in  the  Princeton  oblateness  data 
after  an  estimate  of  static  oblateness  is  removed  (Dicke  1976).  Me  pre- 
fer to  adopt  a somewhat  more  conservative  error  estimate.  The  time  in- 
terval analyzed  with  the  Princeton  oblateness  data  is  97  days  (Dicke 
1974,  1976);  we  will  take  the  uncertainty  in  any  frequency  estimate  to 
be  ±1/2T  (Bendat  and  Piersol  1971).  This  gives  an  expected  error  in 
frequency  of  ±5.15  x 10'3  d"1  or  a range  in  period  of  11.87  to  13.52 
days. 

Me  may  now  estimate  the  probability  that  a peak  with  the  significance 
of  the  12.0715  day  peak  would  occur  by  chance  in  the  interval 
11.87-13.52  days.  This  range  in  period  contains  1,352  spectral  esti- 
mates so  that  we  estimate  the  probability  of  the  12.0715  day  peak  occur- 
ring within  this  range  by  chance  as  (1.352  x 1 0 3 J x (1.4  x 10'6)  or 
-2  x 10'3.  As  we  have  already  indicated,  only  -1/3  of  the  spectral  es- 
timates are  independent,  so  that  calculating  the  probability  estimate  as 
though  all  1352  spectral  estimates  were  independent  produces  a conserva- 
tive probability  estimate. 

The  adjusted  spectral  estimates  in  the  range  11.87-13.52  days  are 
displayed  in  Figure  3 with  confidence  levels  corresponding  to  the  proba- 
bility that  none  of  the  1352  adjusted  spectral  estimates  would  exceed 
the  indicated  values. 

As  a check  on  the  foregoing  analysis,  we  differenced  the  sunspot  data 
and  analyzed  the  differenced  data  in  the  same  manner  as  the  undiffer- 
enced data.  The  analysis  of  the  adjusted  spectral  estimates  produced 
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. similar  results  in  all  three  methods  of  estimating  oz  and  the 

probability  estimates.  Since  the  differenced  data  appears  nearly  uncor- 
l 

related  day  to  day  and  most  of  the  power  at  low  frequencies  is  removed, 
we  performed  an  F test  on  the  least-square  fit  to  the  differenced  data 
at  12.0715  days,  and  obtained  a (chance)  probabilty  estimate  of  6 x 10'7 
(about  a factor  of  2 lower  than  that  from  the  analysis  of  the  cumulative 
distribution  of  the  adjusted  spectral  estimates). 

Me  have  also  split  the  data  into  halves  and  separately  analyzed  the 
first  and  second  half  of  both  the  raw  and  differenced  data.  The  analy- 
sis of  the  adjusted  spectral  estimates  gives  nearly  the  same  results  for 
all  three  methods  of  estimating  oz  for  both  raw  and  differenced  data  for 
the  entire  data  run  and  each  half  considered  separately.  The  products 
of  the  (chance)  probability  estimates  for  the  peak  at  12.0715  days  for 
each  half  considered  separately  are  approximately  equal  to,  but  slightly 
smaller  than,  the  probability  estimates  for  the  entire  data  run  for  both 
raw  and  differenced  data. 
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III.  Discussion 


The  peak  at  12.0715  days,  with  approximately  equal  adjusted  spectral 
estimates  in  the  first  and  second  halves  of  the  data,  is  suggestive  of 
the  influence  of  a stable,  long-lived  periodic  process,  such  as  the  ro- 
tation  of  the  sun's  core.  However,  the  data  could  be  reconciled  with  a 
synodic  rotation  period  of  the  core  which  is  a multiple  of  the  12-day 
period.  In  particular,  a core  rotating  with  a 24  day  synodic  period 
could  cause  an  apparent  12-day  periodicity  if  the  disturbance  produced 
by  the  core  has  not  only  a possible  m = 1 component  but  also  an  m = 2 
component,  where  m is  the  azimuthal  mode  number.  On  the  other  hand,  to 
attribute  a period  of  36  days  or  more  to  the  core  seems  unreasonable, 

i 

since  the  convective  zone  would  then  be  subject  to  decelerating  torques 
from  both  the  core  and  the  solar  wind. 

Me  have  also  split  the  data  into  even  and  odd  activity  cycles,  and 
find  that  the  adjusted  spectral  power  estimate  at  12.0715  days  is  ~7 
times  as  large  for  odd  cycles  as  for  even  cycles.  This  suggests  that, 
if  the  peak  at  12.0715  days  is  due  to  core  rotation,  the  coupling  with 
surface  phenomena  is  probably  magnetic.  Unless  the  relic  magnetic  field 
is  confined  to  one  compact  "active  region",  the  coupling  between  the 
core  and  convective  zone  must  then  have  a strong  m = 2 component.*  Fol- 
lowing this  line  of  argument,  the  sunspot  data  alone  seem  to  favor  the 
interpretation  of  the  12-day  peak  as  being  produced  by  a core  rotating 
with  a synodic  period  of  24.14  days  coupling  magnetically  to  the  convec- 


•indeed  the  presence  of  a peak  in  the  spectrum  at  about  13.5  days  indi- 
cates that  even  the  convective  zone  has  some  kind  of  m = 2 structure. 
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tive  zone.  Me  realize,  however,  that  if  Dicke's  (1976)  interpretation 
of  the  Princeton  oblateness  data  is  correct,  the  12-day  periodicity  of 
the  sunspot  data  is  to  be  interpreted  in  terms  of  a 12-day  core  rota- 
tion. 

Among  the  persuasive  arguments  uhich  Dicke  (1976)  advances  in  favor 
of  his  interpretation,  ue  are  particularly  impressed  with  the  following: 
Only  odd  harmonics  of  the  rotation  period  should  appear  in  the  diagonal 
component  of  the  Princeton  "oblateness"  data  when  the  projections  on 
the  plane  on  the  plane  of  the  sky  of  the  rotation  axes  of  the  sun  and  of 
the  earth  are  aligned.  Tor  an  interval  of  time  satisfying  this  condi- 
tion, Dicke  finds  that  his  data  yield  only  odd  harmonics  for  an  assumed 
rotation  period  of  12.22  days  (sidereal)  but  not  for  assumed  periods 
near  24  days. 

Superficially,  our  findings  seem  to  disagree  with  Dicke's:  Dicke 
(1976)  proposes  that  there  are  two  distortions  per  rotation,  uhich  one 
might  expect  to  correspond  to  a peak  in  the  sunspot  spectrum  at  about  6 
days  rather  than  12  days.  However,  Dicke's  data  are  obtained  by  summing 
signals  from  diametrically  opposed  pairs  of  windows:  hence  what  may  in 
reality  be  a single  localized  distortion  would  in  any  event  appear,  from 
the  Princeton  data,  as  a diametrically  opposed  pair  of  distortions. 

To  summarize,  we  find  a prominent  peak  in  the  sunspot  power  spectrum 
with  a period  of  12.0715  days  uhich  is  consistent  with  Dicke's  (1976) 
12.64  day  period  if  (as  ue  think  appropriate),  the  uncertainty  in  the 
period  inferred  from  the  Princeton  data  is  "±0.8  days.  The  probability 
that  a peak  of  th  s spectral  power  would  occur  within  these  limits  by 
chance  is  estimated  to  be  “2  x 10°. 
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Figure  Captions 

Figure  1.  Least-squares  power  spectrum  generated  as  explained  in  the 
text.  The  median.  70th.  90th  and  99th  percentile  estimated  spectral 
powers  are  plotted  for  each  of  64  equally  spaced  frequency  bins.  The 
sun  symbol  corresponds  to  the  peak  at  12.0715  days. 


Figure  2.  The  quantity  (1-0,  the  probability  that  a spectral  estimate 
will  exceed  S.  as  a function  of  adjusted  power,  S.  The  solid  curve  cor- 
responds to  the  empirically  determined  cumulative  distribution,  the  bro- 
ken line  to  the  cumulative  distribution  for  a variate  distributed  as 
. 5x2  with  two  degrees  of  freedom.  The  right  vertical  axis  is  labeled  by 
the  rank  of  the  corresponding  adjusted  spectral  estimate.  The  numbers  1 
through  5 indicate  the  five  largest  adjusted  spectral  estimates. 


Figure  3.  Adjusted  spectral  power  versus  frequency  for  the  frequency  in- 
terval 7.911  x 10" 2 ± 5.15  x 10'3  day"’.  This  frequency  interval  cor- 
responds to  a range  in  period  of  11.87  - 13.52  days  and  contains  1352 
adjusted  spectral  estimates.  The  confidence  levels  correspond  to  the 
probability  that  none  of  1352  spectral  estimates  distributed  as  .5x2 
with  two  degrees  of  freedom  would  exceed  the  indicated  spectral  power. 
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